Analysis of velocity data

ABSTRACT

The invention is concerned with the analysis of velocity data such as sonic data obtained from wells. The invention includes analyzing time versus depth data and generating velocity versus depth data and generating velocity versus depth function. A significant feature is the recognition that translating the observed time versus depth data along the time axis in time-depth domain can be used to obtain a match with a time versus depth curve corresponding a velocity versus depth function. The invention also makes use of visual displays which show degree-of-fit data in parameter spaces. These displays represent an important means for generating velocity versus depth functions.

FIELD OF INVENTION

This invention relates to the analysis of velocity data such as soniclog data obtained from wells. Inter alia, the invention deals with thegeneration of functions that provide a representation of the variationof seismic velocity in the ground with depth.

BACKGROUND ART

Sonic data are a predominant form of velocity data used in thegeneration of velocity versus depth functions.

Sonic data are acquired from wells which have been bored in the groundby measuring the travel time of sonic signals along segments of theborehole over ranges that usually cover most of its length. Other typesof velocity data may be obtained from an analysis of reflection traveltimes. These reflections arise when the seismic waves from an acousticsource, such as an explosion, i.e. seismic pulse travel through theground and are reflected back to the surface (land or water) where theyare recorded and subsequently processed. Yet other types of velocitydata might be obtained from travel time data of seismic waves emitted byacoustic sources in the surface and picked up by receivers in a boreholeor vice versa. In other cases, the sources may be placed in one boreholeand the receivers in another borehole.

The data, whether from a borehole or surface measurements or from acombination of the two types, are then analyzed and calibrated asnecessary. Such data essentially constitute the observed velocity data.The data are then used to obtain for example a representation of thevariation of velocity with depth in the area to which the data pertainor refer. This representation can be expressed in the form of amathematical function--the velocity versus depth function, for example.

On the basis of the velocity function, it is possible to derive certainestimates relating for example to the depths of various layers in thatarea. These depths are then used to construct a structural model ofrelevant horizons in the subsurface which is then used in determiningfor example the likelihood of the presence of hydrocarbon accumulationsin the ground and the potential volume of such accumulations. The depthestimates also provide a basis for drilling locations and programmes.The velocity representation also provides information on the rockproperties which in turn provide geological information that is relevantto exploring for hydrocarbons.

The variation of velocity in the ground may be representedmathematically by a velocity versus depth or velocity versus timefunction. The mathematical representation gives a smooth description ofthe true variations that tend to be rapidly variable, as will beexplained below. In the present invention, an additional representationis given in the form of a numerical function described in Appendix 1.(In the present text, reference to functions will imply a mathematicaland not a numerical function, unless stated otherwise.)

The instantaneous velocity is a basic quantity in the present invention.The instantaneous velocity at a given depth refers to the speed of theseismic wavefront at that particular depth in the direction ofpropagation or the velocity in an infinitesimally small thickness of therock at that depth which, in mathematical notation, may be representedby

    V.sub.ins =dZ/dT

where Z represents depth and T represents time. It may also be taken asthe velocity value which is measured by the sonic device at that depth.Strictly speaking, the measured velocity in this case represents thevelocity across a finite thickness of rock which is spanned by themeasurement device. However, reference to this velocity as instantaneousvelocity is adequate for most practical purposes and for the purpose ofthe present invention. Sometimes, the term instantaneous velocity isused to refer to the velocity across a smaller or a larger interval thanthat measured by the sonic log device. The present invention includesthese cases regardless of the actual nomenclature used. In particular,it includes cases where interval velocities such as might be obtainedfrom seismic reflection surveys are used as data in the functiongeneration process. These velocities refer to much coarser intervalsthan those measured by the sonic device. The generation of velocityfunctions from interval velocity data and from other data spanningcoarse intervals is given in Appendix 2.

Typically, the value of the instantaneous velocity changes in a rapidfluttered manner down the borehole as illustrated schematically inFIG. 1. In fact, the velocity versus depth function, which is amathematical representation of the variation of the instantaneousvelocity with depth, generally gives a smooth representation thatdescribes the variation in a gross general way as illustratedschematically by the smooth curve in FIG. 1. Hence, reference to thevariation of the instantaneous velocity with depth in terms of thevelocity function would normally imply this smooth variation. Anexception to this smooth variation is the numerical function describedin Appendix 1.

The variation of velocity in the ground may also be represented in termsof instantaneous velocity as a function of time, the time beinggenerally the total travel time of the seismic wave from the surface orfrom some other reference level to the point to which the instantaneousvelocity corresponds. The variation of the instantaneous velocity as afunction of total travel time (not depth) is not readily supportable ongeological grounds. The use of this kind of function stems largely fromthe mathematical convenience of its formulation. The present inventionmainly presents velocity as a function of depth. The generation ofvelocity versus time functions is presented in Appendix 3.

The variation of velocity in the ground may also be represented in termsof slowness functions such as the instantaneous slowness versus depthfunction. Slowness is the reciprocal of velocity and instantaneousslowness is the reciprocal of instantaneous velocity. In mathematicalnotation, it may be represented by S_(ins) =dT/dZ. As in the case ofinstantaneous velocity, the term instantaneous slowness will be used inthe present invention to include situations where the slowness refers toan interval that might be smaller than, equal to or larger than theinterval measured by the sonic device.

The generation of slowness functions such as the instantaneous slownessversus depth function is one of the totally novel concepts introduced bythe present invention. It is implied throughout the present text thatall the aspects relating to the velocity function are also applicable tothe slowness function. Further details of the slowness function aregiven in Appendix 4.

The variation of the instantaneous velocity with depth differs widelyfrom one geological unit to another. A geological unit (or what in thetext is simply referred to as unit) is defined here as a rock formationor layer or the like in which the behavior of the instantaneous velocityas represented by the velocity versus depth function (or any otherfunction describing the variation of instantaneous velocity in theground) is substantially the same throughout its lateral extent. At anyone location, several units will usually be present down to the deepestlevel of interest. In certain situations, the behavior of theinstantaneous velocity within a given unit might vary laterally in agradual manner. Strictly speaking, when the variation is such that thefunction parameters applicable to one location can no longer beconsidered as providing adequate description of the variation ofinstantaneous velocity with depth at another location, the unit shouldbe treated as being a different unit. However, in such cases, it isoften simpler for practical purposes to obtain function parameterestimates at specific locations and then to use these estimates toproduce contour maps of parameter values over the area in question.

Currently, the velocity versus depth function (in the sense of themathematical representation of the instantaneous velocity variation withdepth as described above) is the most commonly used form to representthe velocity variations in the ground. For this reason, the variousfeatures of the present invention will be described using this functionas the main example. However, many other types of function can beemployed to represent the velocity variation in the ground. For example,the variation may be represented by a time versus depth function, aninterval velocity versus depth function, an average velocity (where theaverage is measured from the top of the unit or from some other leveland where the averaging may be simple, root mean square (rms) or anotherform) versus depth, slowness (instantaneous, interval or averageslowness) versus depth, among others. These representations are allequivalent to the velocity versus depth function representation. Forexample, the time versus depth function frequently can be deriveddirectly from the velocity versus depth function; the interval velocityfunction is obtained from the velocity over an interval that isconsidered too coarse to give the instantaneous velocity; the averagevelocity function is obtainable from an integration of the instantaneousvelocity; the slowness function is the reciprocal of the velocityfunction and so on. Likewise, the velocity versus time functionrepresentation has equivalent representations in the form of a timeversus depth function, an interval velocity versus time function, anaverage velocity versus time function, a slowness versus time function,among others. The present invention includes all of the above functions.

The velocity versus depth function can be expressed in a variety offorms. Each form is characterized by a set of parameters. For example,the function pertaining to a given unit might be expressed in the linearform

    V.sub.z =V.sub.0 +kZ

where V_(z) is the value of the instantaneous velocity at depth Z. V₀and k are the parameters whose specific values provide a description ofthe variation of the instantaneous velocity with depth in that unit (seefor example J. A. Legge and J. J. Rupnik, 1943, Least squaresdetermination of the velocity function V_(z) =V₀ +k Z, Geophysics, vol.8, p. 356).

Examples of other forms are the Faust relationship

    V.sub.z =A Z.sup.1/n

whose parameters are A and n (L. Y. Faust, 1951, Seismic velocity as afunction of depth and geologic time, Geophysics, vol. 16, p. 192).

In the form

    V.sub.z =V.sub.0 (1+a Z).sup.1/n

the function has three parameters, V₀ , a and n.

The generation of the velocity function in a given form is essentially aprocess of determining the relevant parameter values. In view of thenon-uniqueness aspect which is discussed later the function generationprocess in the present invention will be considered in the sense ofestimating the parameter value rather than a determination in the senseof obtaining an absolutely accurate or unique value.

A given velocity versus depth function may be transformed into a timeversus depth function or relationship. For example, from the linear formshown above the travel time T corresponding to depth Z can be expressedin the following manner:

    T={log.sub.e  (V.sub.0 +kZ)/V.sub.0 !}/k

Reference to velocity versus depth functions in the present text alsoimplies the corresponding time versus depth functions or relationships.Further discussion of function transformation is presented below inconjunction with function generation domains.

More than one function form may provide an adequate representation ofthe variation of instantaneous velocity with depth in a given unit. Innormal application, however, one form is used per unit which may or maynot be similar to the form used in another unit.

The domain in which the function generation is carried out can, forexample, be time-depth, velocity-depth, velocity-time, slowness-time.The choice of the most suitable function generation domain depends onthe details of the function being generated and on the available data.For example, if sonic data are available for generating a velocityversus depth function then the time-depth domain is the most appropriatedomain. If only interval velocity data from a seismic reflection surveyare available then the velocity-depth domain would often be the mostsuitable domain for generating this function.

Clearly, the data used for the function generation (the observed data)need to be the same as the function generation domain. Likewise, thecurve corresponding to the function being generated (the calculatedcurve) must also be transformed to the form required for thatdomain--Obviously no transformation is required if the function isalready in that form. For example, for generating the linear velocityversus depth function

    V.sub.z =V.sub.0 +kZ

(i.e. for determining V₀ and k) in the velocity-depth domain use theabove form directly and in the time-depth domain use

    T={log.sub.e  (V.sub.o +kZ)/V.sub.o !}/k

These are all equivalent forms of the same function. However, not allfunctions can be readily transformed explicitly to the form required ina particular domain, thus limiting the possibility of generating thefunction in that domain. For example, the function

    V.sub.Z =V.sub.0 +kZ.sup.1/n

cannot be transformed to a time-depth form in a readily usable manner.

In the present invention, function generation is presented mainly in thetime-depth domain but also in the velocity-depth and velocity-timedomains. In the last two domains, the term velocity is used in the senseof instantaneous velocity including the case where the instantaneousvelocity is approximated by interval velocity. However, it is impliedthroughout the present text that all the aspects relating to thevelocity-depth and velocity-time domains are also applicable to theslowness-depth and slowness-time domains as well as averagevelocity-depth and average velocity-time domains.

In the present text, function derivation in the time-depth domain willbe described in detail; derivation in the velocity-depth andvelocity-time domains will be given afterwards. The term time versusdepth curve will be used here to denote a curve in the conventionalsense or a set of discrete time-depth data points with or without acurve joining them.

Function generation techniques known to date have been quiteproblematic. Derivation in the velocity-depth (and velocity-time) domainis usually subject to the scatter of sonic velocities as mentionedabove. A general method for derivation in the time-depth domain had notexisted prior to the present invention. In the derivations that had beenattempted, the time versus depth data had remained fixed in thetime-depth space. This is a fundamental error. Generally, therefore, amatch between a given observed time versus depth curve and a time versusdepth curve corresponding to a velocity versus depth function could notbe effected correctly. For the same reason, time versus depth curvespertaining to the same unit but originating from different wells couldnot be fitted to the same function in a simple way. As explained below,the time versus depth data can generally remain fixed only in thetop-most unit, that is the unit immediately below the surface or,equivalently, under the assumption that the geological section from thesurface to the level of interest consists of one and the same unit.

SUMMARY OF THE INVENTION

The present invention is concerned primarily with a method of analyzingtime versus depth data and generating velocity versus depth functions.The present invention eliminates the aforesaid problems in functiongeneration by demonstrating that translating the observed time versusdepth data along the time axis in the time-depth domain is valid and byusing this translation in order to obtain a match with a time versusdepth curve corresponding to a velocity versus depth function. Thetranslation consists essentially of subtracting or adding an appropriateconstant time value to all the points defining the time versus depthcurve.

The present invention is also concerned with analyzing time versus depthdata, velocity versus depth data, velocity versus time data, slownessversus depth data and slowness versus time data and generating velocityversus time functions, slowness versus depth functions and slownessversus time functions.

A main feature of the present invention is the construction of visualdisplays in the parameter space, preferably in the form of contours.These displays can refer to any set of values particularly valuesdescribing the degree-of-fit between two sets of data such as values ofa discrepancy function or measure or an objective function. The contoursor alternative visual displays delimit regions of parameter values inwhich the value of discrepancy or another degree-of-fit measure islimited by the specific value of each contour. This feature allows apowerful means for generating the velocity versus depth function (andother functions) and for controlling the quality of the parameterestimation process as well as a number of other geological applications.In one aspect of the present invention, this feature is claimed as afundamental part of the invention. It includes the following aspects inparticular singly or in combination:

(a) Any goodness-of-fit criterion (such as a discrepancy function ormeasure) in any appropriate form that expresses the degree of fitbetween the curves being compared or fitted. For example, thediscrepancy in the time-depth domain may be given by a weightednormalized sum of squares of residuals along the time or the depth axisbetween the observed and the calculated function curve.

(b) Space defined by any number of parameters. Although the feature isbest suited for a two-parameter space (i.e. two-dimensional displays) itis also applicable to spaces of higher dimensions.

(c) All pertinent procedures regardless of whether the data analysis,data fitting, discrepancy calculations or the like are carried out inthe time-depth domain, the velocity-depth domain, the velocity-timedomain the slowness-depth domain or the slowness-time domain. Thepreferred domain is the time-depth domain.

(d) All procedures where the parameter estimation process uses theprinciples of these displays implicitly without a visual output.

(e) Output of the value of the discrepancy measure or thegoodness-of-fit criterion or the like in any visual form such as whenthis value is shown in accordance with colour variations, colourintensity, shading, contours, grid values and the like. The preferredform is colour-coded contours.

(f) Output in any form particularly on a monitor screen, screen of anykind, and hard copy of any kind.

In the description of the present invention, only grid/contour displaysin a two-parameter space is presented for the sake of simplicity but theabove inclusions are implied throughout the description as appropriate.

The present invention also includes procedures for generating numericalfunctions, detailed in Appendix 1.

The invention will be described now by way of example only withparticular reference to the accompanying drawings. Function derivationin the time-depth domain will be given as the main example.

BRIEF DESCRIPTION OF THE DRAWING

FIG. 1 is a schematic representation of the variation of instantaneousvelocity with depth as explained already.

FIGS. 2a and 2b are schematic diagrams helpful in describing theillustrating general principles of the present invention

FIGS. 3 to 9 and FIGS. 12 to 14 are plots of velocity vs. the reciprocalof time for some illustrative examples of the present invention, and

FIGS. 10 and 11 are plots of velocity vs. the reciprocal of time andsuperimposition for other illustrative examples of the presentinvention.

FIG. 15 is a time versus depth curve for numerical function generation.

DESCRIPTION OF THE PREFERRED EMBODIMENT (A) FUNCTION DERIVATION IN THETIME-DEPTH DOMAIN

Referring to FIG. 2a there is shown schematically a section through theground beneath the surface (10). Two wells (A) and (B) are formedvertically in the ground at the position indicated. As can be seen theground has a sub-structure comprising a low velocity unit (12) and ahigh velocity unit (14). Below that there is a unit identified as unit(S1).

The graph in FIG. 2a shows time versus depth data obtained in the twowells (A) and (B) for the unit (S1). The curve (A) represents the sonicdata obtained from well (A), that is, at each depth point, the curvegives the total traveltime from some reference level (such as the groundsurface) to that point. Curve (B) represents similar sonic data obtainedfrom well (B). As can be seen both curves are for the same unit (S1).The instantaneous velocity in unit (S1) is assumed to vary as a functionof depth in that unit. Hence at each depth point the slope of thetangents to both curves will be equal, that is to say the two curves arein essence parallel. The shift along the time axis between the twocurves (A) and (B) arises from the differences in travel time in thelayer or layers above the unit (S1). That is to say the shift arisesfrom effects not associated with the unit (S1). If one of the curveswere to be shifted laterally, that is, translated along the time axisrelative to the other until they coincide, the function describing thevariation of instantaneous velocity with depth would not be affected.

The value of the instantaneous velocity at each point along either curvewould not change.

It is the translation of the time versus depth curves along the timeaxis which is a principal foundation of the present invention andprovides the basis for fitting and combining time versus depth curvesfrom a number of wells and from other sources of velocity data. In thepresent invention curves that cover different depth ranges aretranslated so as to fit over their overlapping portions along the depthaxis in accordance with an appropriate goodness-of-fit criterion. Thismay for example be a least squares fit. In particular, the translationprocess provides the basis for fitting observed time versus depth curvesto time versus depth functions (and implicitly their correspondingvelocity versus depth functions). This involves translating a curverepresenting the observed data along the time axis to produce a fit witha time versus depth curve corresponding to an appropriate velocityversus depth function. In FIG. 2a either or both of the curves can betranslated along the time axis to fit the time versus depth curvecorresponding to a given velocity versus depth function which is shownby the curved line (20). If the fit is satisfactory then this functioncan be assumed to provide a description of the variation ofinstantaneous velocity with depth in unit (S1).

FIG. 2b shows a schematic example where a time versus depth observedcurve cannot be satisfactorily fitted to a function curve (Function 1for example) in its original stationary position. It requirestranslating along the time axis. to the dotted position in order toproduce a satisfactory fit with a function curve (Function 2 forexample). Function 2 can then be assumed to provide a description of thevariation of instantaneous velocity with depth in the unit to which thecurve pertains.

It should be noted that the translation technique is not required in thetop-most unit shown for example in FIG. 2a. Time shifts between thecurves for that unit cannot be attributed to differences in travel timesabove that unit. Therefore, strictly speaking, time versus depth curvesfrom different wells and other sources should be substantiallycoincident in the top-most unit.

There will now be given an illustration of the concepts behind thevarious aspects of the present invention, including the visual displaysof discrepancy, using a number of examples.

It will be assumed throughout the text below that reference to observedtime versus depth curves, velocity versus depth curves or velocityversus time curves implies curves pertaining to one given unit unlessstated otherwise.

Consider the case of a unit in which the instantaneous velocity variesexactly in accordance with a specific velocity versus depth function andthe parameters of this function are unique. That is to say no other setof parameters defining the same function type would produce a functionthat would exactly fit the time versus depth data corresponding to thatunit.

As a specific example consider a unit such as S1 in FIG. 2a in which theinstantaneous velocity varies exactly linearly with depth and hasparameter values of V₀ =2000 m/sec and k=1.0 sec⁻¹. This unit is assumedto be 500 meters thick with its top at 1000 meters below the surface.The observed time versus depth curve is assumed to be created frommeasurements made at depth intervals of 20 meters. The observed curve isfitted to a whole series of time versus depth curves calculated fromrespective linear velocity versus depth functions, each of thesefunctions corresponding to a given combination of V₀ and k values. Foreach such combination, translation along the time axis is carried out asnecessary, that is to say the observed curve is translated so as toproduce an optimum fit with the function curve it is being comparedwith. A quantity F(V₀, k) gives a measure of the discrepancy between theobserved curve and each of the functions to which it is fitted. Thisquantity F is represented by the following equation: ##EQU1## whereT_(i) is the observed time at the ith depth point after translation andC_(i) is the function time at the ith depth point, m being the number ofthe digitised points on the observed curve. This equation gives adiscrepancy function, F (V₀ , k), representing what we shall refer tohere simply as the average discrepancy of the fit. Using thisdiscrepancy representation it is possible to produce an array of valuesof F for a range of V₀ and k values. These values are shown in FIG. 3 asa grid/contour display. The contours shown in FIG. 3 are overlaid on thegrid values. It can be seen that the only combination of V₀ and k valuesthat produce an exact fit with the observed time versus depth curve arethose corresponding to the true values of 2000 m/sec and 1.0 sec⁻¹.However, it can be seen that the contour of value 0.5 ms (milliseconds)extends over a range of V₀ values between approximately 1400 and 2600m/sec and k values of between 0.5 and 1.5 sec⁻¹. Any point (i.e. anycombination of V₀ and k) inside this contour will produce a fit with theobserved curve to better than 0.5 ms average discrepancy.

Another observable feature on FIG. 3 is the degree of parallelismbetween the contours. This feature arises as a consequence of theparticular formulation of the function being used (linear variation ofvelocity with depth in the present case).

The average discrepancy values shown in FIG. 3 are based on idealconditions which virtually never arise in practice. The contours showthe possible parameter combination that would produce a fit within agiven margin of uncertainty. They illustrate how, even under these idealconditions, substantial non-uniqueness in the solution for the functionparameters can develop within a small margin of uncertainty. Inpractice, many factors contribute to this non-uniqueness, e.g.observational errors, length of the observed curve and samplinginterval, depth of the unit, difference in the gross behavior of theinstantaneous velocity with depth from that represented by the function,etc.

FIG. 4 is used as an illustration. The kind of general error encounteredin practice was simulated by a series of random errors of maximumamplitude of ±1.0 ms superimposed on the theoretical time versus depthcurve used in FIG. 3. FIG. 4 shows the resulting contours of thediscrepancy function, F(V₀, k), over a similar range of V₀ and k valuesas in FIG. 3. As expected, the random errors produce a general increasein the value of discrepancy along the main solution trough, i.e. ageneral increase in the non-uniqueness of the solution.

The above considerations hold for any set of parameters to varyingdegrees, depending on the specific formulation. Consider a theoreticalexample of a unit, 500 m thick with its top at 1000 m subsurface and inwhich the instantaneous velocity varies in accordance with Faustrelationship, A=1000 (metric) and n=6. FIG. 5 shows a grid/contourdisplay of the values of F(A,n), calculated as in the equation used toproduce FIGS. 3 and 4. As in FIG. 3, the only combination of A and nvalues that produces an exact fit are the true values at A=1000 unitsand n=6 while the substantial non-uniqueness, within only 0.5 ms averagediscrepancy tolerance, is quite evident.

FIG. 6 corresponds to the same data as in FIG. 5 but with the curvelength now increased from 500 m to 1200 m. There is a clear shrinkage inthe extent of non-uniqueness as a result of the increase in curvelength.

It should be noted that the non-uniqueness discussed in this section canbe turned to great advantage in function construction, as will bedetailed below.

Examples of Function Generation Options

The present invention encompasses many options for generating velocityfunction parameters on the basis of the aspects presented above. Theoptions listed below represent specific examples. They will be describednow for function fitting in the time-depth domain. In describing theseoptions, the following will be implied: (a) The time-depth data setshave already been divided according to individual units, i.e. thefunction being generated refers to a specific unit. (b) When differentobserved time versus depth curves are translated along the time axisrelative to one another the translation is carried out to an optimum fitposition. (c) Each value of discrepancy between the observed time versusdepth curve (single or combined) and the time versus depth curvecorresponding to a given velocity versus depth function is obtainedafter translating the observed curve along the time axis to an optimumfit position.

Option 1: Contour Overlap (Multi-Display of Contours): This option isbest suited for functions that have two parameters, e.g. V₀ and k in thelinear case, A and n in the case of Faust relationship, etc.

In this option a contour display, similar to those presented above isgenerated individually for each observed time vs depth curve in turn. Toproduce each display, one parameter is varied along one axis while theother parameter is varied along the other axis. At each parametercombination, a discrepancy value between the calculated time vs depthcurve corresponding to that combination and the observed curve isdetermined at the optimum position after translation. Contours ofdiscrepancy values are then generated in this two-parameter space. Thediscrepancy function may be defined in accordance with any appropriatecriterion. The choice of the contour (discrepancy) value that isconsidered to be acceptable for the pertinent observed curve depends onthe details of the problem such as data accuracy, curve length, etc.This contour may be regarded as corresponding to the "maximum tolerance"value for that observed curve. The "maximum tolerance" value need not bethe same for all the observed curves. Each point contained inside thecontour of "maximum tolerance" corresponds to a function that fits theobserved curve to better than the discrepancy value of the contour.

The "maximum tolerance" contours corresponding to the individualobserved time vs depth curves are then superimposed such that theparameter axes coincide at the corresponding values on all displays.Each point in the zone of overlap between the "maximum tolerance"contours will then correspond to a function that fits all the timeversus depth curves to better than the respective tolerance value ofeach curve.

In practice, once the overlap had been obtained, contours ofsuccessively lower discrepancy values could be attempted. This wouldcause the zone of overlap to shrink in size but, as long as the zone hadnot vanished completely, the points within this zone will correspondincreasingly to closer fits with the observed time versus depth curves.Any point within this zone (such as a central point) can then beselected to represent the parameters for that unit.

Conversely, if one or more of the "maximum tolerance" contours do notoverlap, the value of tolerance may be increased to test the level oftolerance required in order to obtain overlap. Such tests can helpdecide on the validity of certain assumptions, e.g. whether the datafrom the various wells/sources truly correspond to the same unit.

FIGS. 7 to 11 serve as an illustration. Each of FIGS. 7, 8, and 9corresponds with an observed time versus depth curve. Each curve hasbeen synthetically generated for a linear velocity function with V₀=2000 m/sec and k=1.0 sec⁻¹ and with random errors of ±1.0 ms maximumamplitude superimposed on the time values. A sampling interval of 20 mis used in all cases. The curves are designated Curves 1, 2 and 3respectively. The other details are as follows: Curve 1: Depth range800-1200 m. Shifts applied to the times, decreasing progressively from0.0 ms at 800 m to -10 ms at 1200 m (to simulate shape distortion, i.e.a departure from a truly linear function). Curve 2: Depth range1300-1900 m. Shifts applied to the times, increasing progressively from0.0 ms at 1300 m to +12 ms at 1900 m (producing a distortion in theopposite direction to that of Curve 1)

Curve 3: Depth range 2000-2500 m. No shape modification.

The 1 ms, 2 ms and 3 ms discrepancy contours are shown on all threedisplays (FIGS. 7, 8 and 9). A tolerance level of 2 ms producesconsiderable overlap between the contours of Curves 1 and 3. Reducingthe tolerance to 1ms still preserves significant overlap. FIG. 10 showsa superimposition of the 1 ms contours of the three displays. The areaof overlap between Curves 1 and 3 (dotted) delineates the region in theV₀, k space where every point would produce a function that fits bothcurves to better than 1 ms average discrepancy. If only these two curveswere being considered then the parameters for the unit could be selectedat the point V₀ =2500 m/sec and k=0.75 sec⁻¹, for example.

The 1 ms contour corresponding to Curve 2 does not overlap with theother two. In fact, a 2 ms tolerance level would still not produce anoverlap. FIG. 11 shows a superimposition of the 3 ms contourscorresponding to Curves 2 and 3 and the 2 ms contour corresponding toCurve 1. At this level of tolerance, there is overlap. If it isconsidered that this level of tolerance is acceptable for the particulardata sets in hand then any point within the overlap area (dotted) wouldgive suitable parameters for the unit. For example, the point V₀ =2500m/sec, k=0.65 sec⁻¹ may be selected as representing the parameters ofthe unit.

Option 2: Single Display of Contours: As in the previous option, thisoption is best suited for a function that depends on two parameters.

One way of producing a single display of contours could be achieved forexample as follows: All the time versus depth curves are translatedalong the time axis until they fit optimally along their overlappingrange in depth. If one curve does not overlap with the other(s) or ifthe overlap is not considered to be adequate then an appropriatefunction needs to be temporarily generated for that curve andextrapolated to produce the required overlap.

Obviously, if the entire data set consists of one curve only then notranslation is applicable at this stage. Also, no translation isrequired if the observed data pertain to the top-most unit.

When the translation process is completed, the combined curves can betreated as a single curve. This curve will generally be multi-valued intime at every digitised depth point over the range of overlap. The curvewill cover the same or greater depth range than any of the constituentcurves. A contour display corresponding to this single composite curveis then produced, in a similar manner to those described above. Thediscrepancy function on which such display is based may for example bein the following form ##EQU2## where r is the number of observed curves,m(j) is the number of points on the jth curve, N is the total number ofall observed points being compared, T_(ij) is the observed time, aftertranslation, on the ith depth point of the jth curve, C_(i)(j) is thefunction time at the ith depth point of the jth curve and W_(j) is aweight assigned to the jth curve.

The parameter values are then chosen from within the region delineatedby the contour of lowest discrepancy such as the parameter values at thepoint producing the lowest discrepancy value.

FIG. 12 is a grid/contour display based on a combination of the curvescorresponding to FIGS. 4 and FIG. 8 (Curve 2). The two curves overlapover the depth interval of 1300 m and 1500 m. Both curves are weightedequally at 1.0. The total length of the combined curve is 900 m.

The increased length of the time-depth curve leads to a reduction innon-uniqueness (c.f. FIGS. 4 and 8). The discrepancy values along themain trough are higher than those in FIGS. 4 and 8. This is due to thedistortion on Curve 2 which, unlike the situation on FIG. 8, must nowalso be accommodated by the somewhat different curve that corresponds toFIG. 4. Despite the distortion on Curve 2 and the imposed random errors,the discrepancy along the main trough is still reasonably low. For achoice of parameters of V₀ =2200 m/sec and k=0.85 sec⁻¹, a fit with anaverage discrepancy of less than 2.0 ms is obtainable.

Another way of producing a single display of contours could be achievedfor example without combining the curves. Each time versus depth curveis treated individually. For each curve one parameter is varied alongone axis while the other parameter is varied along the other axis and adiscrepancy value is worked out at each parameter combination. Thediscrepancy function may be of any appropriate form. The same range ofparameter values is scanned for every time versus depth curve. For eachparameter combination a value for the overall discrepancy is thencalculated. This value could be represented for example by ##EQU3##where G_(j) is the discrepancy value corresponding to the jth timeversus depth curve, W_(j) is the weight assigned to the jth curve and ris number of observed time versus depth curves. A contour display ofoverall discrepancy values is then produced in the two-parameter space.

There is a possible third procedure that draws from the above twoalternatives. This procedure consists of combining some of the curves toform one composite curve, other curves to form a second composite curveand so on. Each of the composite curves (and any curves that remainuncombined) can then be treated individually as a single curve. Theoverall discrepancy is then worked out at all parameter combinations anda single display of contours is then produced.

The three alternatives provide a useful choice, each alternative beingmore suited to a given field and geological setting than the other two.

Option 3: Direct Parameter Determination: This option can be used withany function, regardless of the number of parameters on which itdepends. It is particularly suited to functions depending on more thantwo parameters.

One way of achieving direct parameter determination consists oftranslating the observed time versus depth curves along the time axisuntil they fit optimally between themselves to obtain a combined curve,in the manner described in Option 2 for the case involving curvetranslation. The parameters are then worked out direc tly from thecomposite curve without producing a visual display. The process normallyinvolves the search for a minimum of an appropriate objective function.A suitable objective function would be in the form of a discrepancyfunction similar to that presented in Option 2 for the case of curvecombination. The search process would usually be carried out throughsome form of optimization technique. The parameter combinationcorresponding to the minimum objective function would then be regardedas the required solution.

An alternative to the non-linear search is to scan a range of parametervalues in a similar way to that used in the previous two options (butwithout producing a visual display). It can then select the parametercombination producing the least discrepancy value or could delimit theregion in the parameter (hyper)space where the discrepancy values arewithin a specified discrepancy, for example. This alternative issuitable when the optimum parameter values are known fairlyapproximately. However, it becomes impractical when the number ofparameters exceeds 3 in view of the large number of discrepancyevaluations that is needed.

The direct parameter determination process can also be carried outwithout combining the curves. Each of the observed time versus depthcurves is treated individually. At each iteration during theoptimization process, a time versus depth curve corresponding to theparameter combination at that iteration is generated (i.e. calculated).Each observed curve is then translated along the time axis to a best fitposition with the calculated curve. The discrepancy between thecalculated and each of the observed curves is then added in a similarmanner to that used in the corresponding case of Option 2. The overalldiscrepancy would thus constitute the objective function. Theoptimization process continues in this way until convergence to anoptimum has been achieved.

(B) FUNCTION DERIVATION IN THE VELOCITY-DEPTH DOMAIN

The general aspects discussed in connection with the time-depth domainapply to the velocity-depth domain. For example, let us take the case ofthe linear velocity versus depth function

    V.sub.z =V.sub.0 +kZ

the main symbols being as described earlier. Because of its particularformulation, this function corresponds to a straight line which has aslope equal to k and an intercept at Z=0 equal to V₀. FIG. 13 is basedon the same data as FIG. 3 but is constructed in the velocity-depthdomain. That is to say, the observed velocity-depth data pertain to anideal unit in which the instantaneous velocity varies exactly linearlywith depth with V₀ =2000 m/sec and k=1.0 sec⁻¹, the unit being 500 mthick with its top 1000 m deep. FIG. 13 shows a grid/contour display ofthe average discrepancy between the input (observed) velocities and thestraight line that corresponds to each particular combination of V₀ andk at each grid point. Non-uniqueness within a comparatively small levelof tolerance is clearly displayed. Note also the exact symmetry in allstraight line directions with respect to the exact solution at V₀ =2000m/sec and k=1.0 sec⁻¹, due to the particular form of this function.

Under the ideal conditions of FIG. 13 it may seem that working in thevelocity-depth domain should be at least as desirable as working in thetime-depth domain. However, in practice, working in the velocity-depthdomain is fraught with difficulties due to the rapidly fluctuatinginstantaneous velocity values in the familiar manner of a sonic log.FIG. 14 was produced, in the velocity-depth domain, after superimposingrandom errors of maximum amplitude of ±1.0 ms on the time data used inFIG. 13. As FIG. 14 shows, the addition of errors here has greatlyincreased the discrepancy level (note that the units of discrepancy are10 times greater than those of FIG. 13). The addition of errors has alsoconsiderably flattened the solution trough. Clearly then, the derivationof parameters from this kind of data in this domain is quite unreliable.However, working in the velocity-depth domain is preferable to workingin the time-depth domain in certain cases, such as when the availabledata are already in velocity-depth format and in which the velocity doesnot show significant flutter.

In the velocity-depth domain the discrepancy between the observedvelocity versus depth data and the velocity versus depth function curvecould again be given by any appropriate measure such as a weightednormalized sum of squares of residuals along the velocity or the depthaxis. Unlike the case in the time-depth domain, in this domain, thegeneration of the velocity versus depth function does not require thetranslation of curves either between themselves or onto the functioncurve.

Otherwise, the application of all the options described in conjunctionwith the time-depth domain follows in a straightforward manner. Forexample, suppose that sonic data for a specific unit are available froma number of wells in an area. Suppose that an individual contour display(similar to that of FIG. 13) corresponding to each well is produced inthe velocity-depth domain. In this case, the contour overlap procedure(Option 1) can be applied in the same way as that used in the time-depthdomain. If, instead of treating the wells individually, the velocityversus depth data from all the wells are considered as part of the samedata set then a single contour display (Option 2) can be produced. Thefunction can also be generated simply in Option 3 through a linearregression fit of the observed data onto the function curve (a straightline in this case); the function parameters in this case being givendirectly by the intercept and slope of the straight line.

(C) FUNCTION DERIVATION IN THE VELOCITY-TIME DOMAIN

Function generation in this domain can be useful for example when theobserved data are already in the form of velocity versus time data suchthat might be obtained from a seismic reflection survey.

If the observed data are not in the form of velocity versus time datathen they should be transformed to this form. Similarly, the functioncurve should be calculated in this form. The travel time is measuredeither from the surface or from some other reference datum.

As in the case of the time-depth domain, for a unit below the top-mostunit, any time shift between the different observed curves isattributable to differences in travel time outside the unit itself.Thus, the observed curves could be translated along the time axisrelative to one another and relative to the function curve.

The generation of the velocity versus depth function proceeds in thisdomain in the same way as in the case of the time-depth domain: Thetranslation of the observed curves along the time axis to fit optimallybetween themselves (if curve combination is required) and with thefunction curve, is carried out as appropriate. Any suitable measure ofdiscrepancy between the observed data and the calculated function curvemay be used. Options 1, 2 and 3 can be used in the same way as describedin the case of the time-depth domain.

MISCELLANEOUS COMMENTS

The choice of the function form to be adopted for a given area isprimarily a question of how well would such function fit the data. Inthe case where only one time vs depth data set is available, the maincriterion is the magnitude of discrepancy produced by the best fittingfunction, in that function form. For example, suppose that a time vsdepth curve for a particular unit in a particular area can be fitted toa specific Faust function to within a discrepancy of 1 ms. Suppose thatthe curve could not be fitted to other tested functions including thelinear velocity function to better than 3 ms. Clearly, in this case, theFaust relationship would be the more suitable function form. In the caseof more than one data set, other considerations come into play, e.g. ifusing Option 1, the function form ought to be that which produces acontour overlap at a reasonable discrepancy level that cannot beimproved upon with the use of another function form; if using Option 2,the function form should produce a better overall discrepancy than ifanother function form was used, etc.

Option 1 is an extremely efficient and versatile option. It can also beuseful in indicating the possibility of lateral changes in the unit.Option 2 works particularly well in areas where the unit does not suffersignificant lateral changes. Both Options 1 and 2 are highly suited tointeractive workstation and PC applications. They are also suited tonon-interactive environment, especially Option 2. Option 3 does notenjoy the overview and control of Options 1 and 2. However, it can bereadily automated for non-interactive use. It is also the obvious choicewhen the formulation involves more than two parameters.

In the parameter determination process the above or similar options aswell as the domains used may be combined as appropriate, providing thatthe parameters being estimated refer to the same function form. Forexample, consider a situation where seismic reflection data and sonicdata from two wells are the only available data in a given area. Supposethat the reflection data are used in the velocity-depth domain toproduce a contour display via the single contour display option. Supposealso that the sonic data are used in the contour overlap option via thetime-depth domain. In this case the contour display from the reflectiondata could be used in conjunction with the sonic data for an overallparameter determination through the contour overlap option.

The present invention includes all situations where the time-depth,velocity-depth, velocity-time, slowness-depth or slowness-time data usedin generating the function are subjected to various corrections to allowfor such factors as dip, structural configuration, changes in the depthof burial (reduced overburden), overpressuring, datum changes,lithological changes, etc. The various procedures given in the presentinvention can be adapted, if need be, to suit the requirements of eachof these problems. As an example, consider for a given unit a time-depthdata set from a certain well. Suppose that the unit at the well locationis known to have undergone a reduction in overburden of, say, 800 m,with respect to its maximum past depth of burial. In this case, 800 mcan be added to each depth value in the data set (with a possiblefurther adjustment for decompaction). The corrected time-depth data setcan then be used in function generation. If the amount of overburdenreduction is not known then an application of the multi contour displayfor example can help to estimate that amount. In fact, the proceduresdescribed in the present invention provide the possibility to estimatethe relative uplift or relative overburden changes for a given unitbetween any two locations for which velocity data are available.

Another example of adapting the above procedures relates to depthreference level or datum. In the description of the present invention,depth was referred to the surface. However, the presented methodologiesapply equally to cases where the depth is referred to some other levelsuch as sea level, sea bed, or the top of a particular unit. The onlyadaptation that is required in such cases is to adjust the depths bysubtracting or adding an appropriate value. (Although in the majority ofcases in practice, the depth should be referred to the ground surface,there are certain applications where reference to some other levelsimplifies the treatment whilst still remaining within the accuracyrequirements of the problem.)

APPENDIX 1 GENERATION OF NUMERICAL TIME VERSUS DEPTH AND VELOCITY VERSUSDEPTH FUNCTIONS

For the present purpose, a numerical function means one that is definedby a series of discrete time versus depth, velocity versus depth orvelocity versus time points. The velocity is defined as being dependenton depth (but not on travel time below the top-most unit) but noassumption is made as to whether the variation of velocity with depthfollows any systematic or unsystematic pattern--The rare case wherevelocity is defined as being dependent on travel time (but not on depth)is briefly covered in discussing function generation in thevelocity-time domain. As in the case of analytical functions, the termvelocity will be used in the case of numerical functions mainly in thesense of instantaneous velocity but will also encompass cases where theinterval to which the velocity refers is relatively coarse, i.e. thevelocity is interval rather than instantaneous velocity. The equivalenceof average velocity and slowness (instantaneous, interval and average )will also be implied. It will be assumed that the input observed datapertain to the same unit. The term curve will be used to denote a curvein the conventional sense or a set of data points with or without acurve joining them. The preferred form of the function is a time versusdepth function.

Function Generation in the Time-Depth Domain

The procedure works on time versus depth data. Well sonic data shouldpreferably be sampled at sufficiently frequent intervals to provideadequate representation of the variation of the velocity in the ground(e.g. at 1 foot intervals). Data from seismic reflection sources(usually given in velocity versus time values) should be suitablytransformed before input. The generated function is then given by theset of time versus depth values that best fit the translated curves inaccordance with a suitable goodness-of-fit criterion. The following isan example of how the generation of the function might be carried out:

1. If the data pertain to the top-most unit, the time versus depthcurves from various wells and other sources (the observed curves) shouldbe fairly coincident. No curve translation is required in this unit. Thenumerical function may then be generated in the form of a time versusdepth data set which best fits the translated data for example on thebasis of a weighted least squares criterion.

If all the observed curves are sampled at the same depths then the timeversus depth function curve should preferably be sampled at those samedepth values. If the observed curves are not sampled at the same depthsthen some form of interpolation, such as a linear interpolation, will beneeded.

If the observed time versus depth curves are significantly divergentthen that may indicate lateral variation in the unit or largemeasurement errors. In such cases, a numerical function for thistop-most unit cannot be reliably generated.

2. In the case of deeper units than the top-most unit, the observedcurves are translated along the time axis relative to one another untilan optimum fit is obtained between the translated curves. Thetranslation may be carried out in any order such as in order of curvelength. The goodness-of-fit criterion could be for example least squaresalong the time axis. FIG. 15 shows schematically three observed curvesin the translated position.

3. If the curves do not overlap adequately (e.g. have only an overlapover a depth range equal to 10% of the length of the shorter curve) ordo not overlap at all then an analytical function curve (e.g. onecorresponding to Faust relationship) is used. This curve is fitted toeither the stationary curve (or curves) or to the curve being translatedor to both. The fitted curve is then extrapolated according to thisfunction so as to produce adequate overlap (e.g. 50% of the length ofthe shorter curve).

4. The numerical function may then be generated in the form of a timeversus depth data set which best fits the translated data on the basisof, for example, a weighted least squares criterion. The sameconsiderations regarding the depths at which the observed curves aresampled apply as in the case of the top-most layer. The generatedfunction is shown schematically in FIG. 15.

It should be noted that abrupt breaks can occur in the function at theextremity points of each curve as the curve starts to be included in theaveraging process at its shallow extremity or finishes from beingincluded in the averaging process at its deep extremity. This is shownschematically in FIG. 1 i5=. These kinks may be smoothed out through theapplication of a suitable spline fit or by distributing the magnitude ofthe breaks both upwards and downwards in a gradual manner. In all cases,the smoothing should ensure that successively deeper points haveincreasing time values in order to avoid the implied negative velocityotherwise.

Function Generation in the Velocity-Depth Domain

The procedure works on velocity versus depth data. No curve translationis applicable in this domain. The generated function is given by a setof velocity versus depth values that best fit the observed curves inaccordance with a suitable goodness-of-fit criterion. The following isan example of how the function generation might be carried out:

1. At each depth point, a weighted average of all the velocity valuesoccurring at that depth is calculated--The number of velocity valuesoccurring at each depth point is equal to the number of curvesstraddling that depth point. The average value becomes the functionvalue at that depth. As in the case of function generation in thetime-depth domain, interpolation may be required if the observed curvesare not all sampled at the same depth values.

2. The generated function is then represented by the series ofvelocity-depth points that were derived in the previous step.

3. If there are depth ranges that are not covered by any curve then thegap can be bridged up for example by joining the deepest function pointon the shallower segment with the shallowest function point on thedeeper segment. The two points could be joined with a straight line orsome other form of line as judged to be appropriate.

Function Generation in the Velocity-Time Domain

The procedure works on velocity versus time data such as data whichmight be obtained from seismic reflection sources. The generatedfunction is given by a set of velocity versus time values that best fitthe observed curves in accordance with a suitable goodness-of-fitcriterion. If the velocity is defined as being dependent on depth (butnot on travel time) then the following is an example of how the functiongeneration might be carried out in this domain:

1. As in the case of the time-depth domain, no curve translation isapplicable in the top-most unit. Below that unit, time shifts betweenthe different observed curves are attributable to differences in traveltime outside the unit for which the function is being generated. Thesecurves are then translated along the time axis relative to one anotheruntil an optimum fit between them is achieved.

2. At each time point, a weighted average of all the velocity valuesoccurring at that time is calculated--The number of velocity valuesoccurring at each time point is equal to the number of curves straddlingthat time point. The average value becomes the function value at thattime. As in the case of function generation in the time-depth domain,interpolation may be required if the observed curves are not all sampledat the same depth values.

3. The generated function is then represented by a series ofvelocity-time points that were derived in the previous step.

4. Gaps (i.e. time ranges not covered by any curve) can be bridged up bya suitable form of line joining the deepest function on the shallowsegment and the shallowest function point on the deep segment.

In the rare case where the velocity is defined as being dependent ontravel time (but not on depth) no curve translation will be applicablein any of the units. Otherwise, the function is generated as in steps 2to 4 above.

APPENDIX 2 VELOCITY FUNCTION GENERATION FROM DATA OBTAINED OVER COARSEINTERVALS (REFLECTION AND WELL SHOOT DATA AS EXAMPLES)

Seismic reflection data constitute the main source of interval velocityinformation in seismic work. The basic quantity is that obtained fromthe maximum coherency between the seismic reflections when the seismictraces are stacked. This quantity is known in the geophysical literatureas stacking velocity. The stacking velocity data need to be correctedfor various factors in order that the data should provide realisticestimates of the velocity in various intervals in the ground uponfurther analysis. Interval velocity data may also be obtained throughsuch geophysical techniques as modeling and imaging on 2D and 3D datasets. Generally, velocities from such processes require fewercorrections than stacking velocity data. In stacking velocity work, aunit typically consists of between one and six intervals. Each intervalis typically between 100 m and 700 m in thickness.

Well calibration surveys and vertical seismic profiling techniques areexamples of techniques that produce direct time versus depth data fromwhich interval velocities can also be derived. These techniques usuallywork by letting off acoustic energy from a source (such as dynamite oran air gun) at the surface and recording the signal and travel time at asuccession of levels down a well using one or morereceiver(s)--Sometimes the source and receiver positions are reversed.Each depth level can be regarded as delineating the base of the intervalabove and the top of the interval beneath that level. The data obtainedfrom such techniques will be generically referred to as well shoot data.

The data obtained over coarse intervals will be collectively referred toas coarse interval data. Function generation will refer to thegeneration of velocity (or implied slowness) versus depth functions.Velocity versus time functions are dealt with in Appendix 3.

Function Generation in the Time-Depth Domain

Interval velocities from reflection data are usually available in theform of velocity-time data pairs. Time-depth data pairs can be generatedfrom this data set using basic principles, the depth to the base of themth interval, D_(m), being given by ##EQU4## where V_(i) is the intervalvelocity in the ith interval and T_(i-1) and T_(i) are the travel timesto the top and base of the ith interval.

Well shoot data are usually directly available as time versus depthdata.

For a given unit, each time versus depth data set, whether from surfacereflection data, from a well shoot or from other sources, can be treatedin the same way as time versus depth data originating from sonic logswhich are usually much more finely sampled. The coarse interval data cantherefore be used in any appropriate combination, with or without soniclogs, to generate the function through, for example, Options 1, 2 or 3as described in the text. However, a preferred procedure in the case ofa two-parameter function (e.g. V₀ and k in the linear case), is asfollows:

Use the coarse interval data to produce a single contour display (Option2) based on the alternative whereby the observed curves are comparedindividually (not combined) with the calculated function curve. Thisdisplay can then be used in a contour overlap approach (Option 1) inconjunction with finely sampled data or with other single contourdisplays from determinations in the velocity-depth domain for example.If there are no finely sampled data or other single contour displaysthen the function generation process is regarded as complete; thefunction parameters can be selected from the display as appropriate.

Function Generation in the Velocity-Depth Domain

This domain is particularly suitable for generating functions that havea simpler velocity-depth form than a time-depth form.

For interval velocities from reflection data it is preferable to assignthe interval velocity value of each interval to the mid point of thatinterval and to treat it as the instantaneous velocity value at thecorresponding depth. The depth to the mid point of the mth interval,D_(mid), being given by ##EQU5## where V_(m) and V_(i) are the intervalvelocities in the mth and ith intervals, T_(m-1) and T_(m) the traveltimes to the top and base of the mth interval and T_(i-1) and T_(i) thetravel times to the top and base of the ith interval.

For well shoot data, the velocity in each interval is obtained from theinterval thickness (difference between the depths below and above theinterval) divided by the transit time (difference between the traveltime to the base and top of the interval). Again, it is preferable toassign this velocity value to the mid point of the respective intervaland to treat it as the instantaneous velocity value at the correspondingdepth.

For a given unit, each velocity versus depth data set, whether fromsurface reflection data, from a well shoot or from other sources, can betreated in the same way as velocity versus depth data originating fromsonic logs. The coarse interval data can therefore be used in anyappropriate combination, with or without sonic logs, to derive thefunction through, for example, Options 1, 2 or 3 as described in themain text.

Function Generation in the Velocity-Time Domain

Working in this domain is useful for interval velocity data fromreflection surveys which are usually available in a velocity-time formand do not require transforming.

For a unit below the top-most unit, any time shift between the differentobserved curves is attributable to differences in travel time outsidethe unit itself. Thus, the observed curves can be translated along thetime axis relative to one another and relative to the function curve asnecessary.

The function generation process can follow the same procedure describedin the main text under generation in the velocity-time domain. Thepreferred procedure is that followed in the case of time-depth domaingeneration described in this Appendix.

APPENDIX 3 THE GENERATION OF VELOCITY VERSUS TIME FUNCTIONS

Instantaneous velocity as a function of travel time is difficult tojustify geologically (apart from the top-most unit where the absolutetravel time is related to depth through the instantaneous velocity). Theuse of velocity versus time functions in geophysics stems mainly fromthe mathematical convenience of the formulation since the time-depthrelationship is generally simpler to derive from this type of functionthan from the velocity versus depth function.

In the generation of velocity versus time functions no translation ofcurves along the time axis is applicable in any of the generationdomains. This is because such a translation would invalidate thegeneration process as the velocity is defined as being a function oftravel time. Apart from this constraint, the generation process followsthe same general methodology presented in the main text for the case ofvelocity versus depth functions. Any suitable measure of discrepancybetween the observed data and the calculated function curve may be usedin any of the domains. Options 1, 2 and 3 can be used as required. Thegeneration of the linear velocity versus time function in thevelocity-time domain can be used as an example. This function has theform

    V.sub.T =W.sub.0 +aT

where V_(T) is the instantaneous velocity at time T and W₀ and a are thefunction parameters. The function corresponds to a straight line thathas an intercept at T=0 equal to W₀ and a slope equal to a. A contourdisplay of discrepancies can be produced over scanned ranges of W₀ and avalues. Options 1 and 2 can be used as required. Generating the functionthrough Option 3 can be carried out by means of a linear regression fitof the observed data onto the function curve (a straight line in thiscase), the function parameters being given directly by the intercept andslope of the regression function curve.

APPENDIX 4 THE SLOWNESS FUNCTIONS AND DOMAINS

The term slowness is used in seismic work to denote the reciprocal ofvelocity. It is not a widely used term. To date, slowness has been usedin problems involving certain inverse problems such as seismictomographic inversion. The use of slowness versus depth functions toprovide possible alternative to velocity versus depth functions is atotally novel concept introduced in the present invention. In thepresent Appendix, the term slowness will be used mainly in the sense ofinstantaneous slowness but will also include interval slowness whereappropriate. Also, the term slowness function will be used in the senseof slowness versus depth function unless stated otherwise.

The scope of application of slowness functions is exactly the same asthose of velocity functions, that is, mainly for converting seismictimes to depths and for various geological investigations like reducedoverburden and lithological or facies changes. Indeed, all of theaspects described in connection with velocity functions apply toslowness functions in the same way apart from the fact that velocity andslowness are the reciprocal of one another. The reason for introducingand highlighting slowness in the present invention is that itsignificantly widens the range of possible functions that are availableto the geophysicist. Slowness as a function of depth is fullysupportable geologically on the same basis as the velocity versus depthfunction. This aspect, combined with the ease of mathematicalmanipulation of slowness versus depth functions, makes it possible toeffect the above applications very accurately and conveniently throughthe use of quite versatile slowness versus depth functions.

Examples of slowness functions and the corresponding time-depthrelationships are given below (Z and T represent depth and timerespectively. ##EQU6## is instantaneous slowness. The remaining symbolsare function parameters): ##EQU7##

The above functions are quite versatile and useful while theircorresponding time-depth relationships are simple in form and easy toderive. Thus, functions forms that are quite intractable for velocityversus depth functions can be used in the slowness approach, making itpossible to simulate a wide range of ways in which the velocity variesin the ground.

The generation of the slowness function in the time-depth domaininvolves producing the function curves through the time-depthrelationship corresponding to the function being generated. This step isidentical to that in the case of velocity versus depth functiongeneration. The remainder of the generation methodology, including theuse of Options 1, 2 and 3 is also identical.

The procedure for generating the slowness function in the other domains(velocity-depth, velocity-time, slowness-depth, slowness-time) differsfrom that for generating the velocity versus depth function only in thetrivial transformation (taking the reciprocal) if the slowness versusdepth function is generated in a velocity domain or the velocity versusdepth function is generated in a slowness domain. Otherwise, themethodology is again identical to that in generating the velocity versusdepth function.

For example, let us consider the linear slowness function shown above.In the slowness-depth domain the function curve is a straight line thathas an intercept at Z=0 equal to S₀ and a slope equal to -g . Thediscrepancy between this line and the observed curve or curves for theunit for which the function is being generated can be measured by anyappropriate goodness-of-fit criterion. A contour display ofdiscrepancies can be produced over the scanned ranges of S₀ and gvalues. Options 1 and 2 can be used as required. The contour overlapoption can have a combination of displays generated through differentdomains. Generating the function through Option 3 can be carried out bymeans of a linear regression fit of the observed data onto the functionstraight line. The function parameters are then given directly by theintercept and slope of the regression line.

Generation of the same function in the velocity-depth domain for examplewould require function transformation to

    V.sub.Z = S.sub.0 -gZ!.sup.-1

and the function curve will no longer be a straight line.

The requirement for generating velocity functions in slowness domains orslowness functions in velocity domains is likely to arise in specialcases where the particular function formulation renders it more suitablefor generation that way.

Slowness versus time functions have no analytical advantages over otherfunctions nor are they geologically justifiable. Their use is likely tobe only for special applications where the ground down to the horizon ofinterest is treated as being one unit. The generation of slowness versustime functions follows the same methodology as in generating velocityversus time functions with function transformation being carried out tobe compatible with the domain in which the function is being generated.

I claim:
 1. A method of determining, from well sonic data and/or seismicdata resulting from sonic/seismic pulses derived at one or morelocations, the parameters of a mathematical function representing thevariation with depth of a velocity-related quantity associated withproperties of the ground below the surface to thereby provideinformation relating to geological or related aspects of groundstructure, the method comprising responding to the data to constructrepresentations of degree of fit of said data in parameter spacerelating to the mathematical function.
 2. The method of claim 1 wherein,in the case of time-depth data below a topmost unit for deriving thedata, the representations of degree of fit data are constructed bytranslating the data along the time axis.
 3. The method of claim 1wherein a quantity indicative of (a) said data or (b) output informationis represented as a function of instantaneous velocity or slowness (i.e.the inverse of velocity) or interval or average velocity or slowness ofsonic/seismic vibrations propagating in the ground.
 4. The method ofclaim 3 wherein the quantity is the instantaneous velocity ofsonic/seismic vibrations propagating in the ground at a particulardepth.
 5. The method of claim 3 wherein the quantity is theinstantaneous slowness of sonic/seismic vibrations propagating in theground at a particular depth.
 6. The method of claim 1 wherein therepresentations are displayed.
 7. The method of claim 6 wherein thedisplayed representations include contours of equal degree of fitmeasures.
 8. The method of claim 6 further including the step ofconstructing a plurality of representations of degree of fit data inparameter space, each of the plural representations pertaining to one ormore sets of well and/or seismic velocity data.
 9. The method of claim 1further including determining a value or a range of values for eachparameter.
 10. The method of claim 1 further including displaying withinthe parameter space contours or parameter values within the contours,and using the contours or parameter values for or in a parameterdetermination process.
 11. The method of claim 10 further includingusing the contours or parameter values within the contours to obtaingeological information about the ground wherein the sonic vibrationspropagated in the ground.
 12. The method of claim 1 wherein themathematical function has two or three parameters.
 13. The method ofclaim 12 wherein the mathematical function is the linear function V_(z)=V₀ +kZ or a transformation thereof.
 14. The method of claim 12 whereinthe mathematical function is the Faust function (V_(z) =A Z^(1/n)) or atransformation thereof.
 15. A method of approximately determining firstand second parameters of a theoretical equation having two parametersand two variables describing how an acoustic wave travels as a functionof depth through a unit of the earth, the acoustic wave travellingthrough the unit to a sonic/seismic transducer, comprising responding tosonic/seismic data derived by the transducer responding to the acousticwave by establishing at least one two dimensional contour, each contourrepresenting a value of acoustic wave travel discrepancy from thetheoretical equation of the acoustic wave travel time to the transducerthrough the first unit, presenting the at least one contour on a twocoordinate axis display having a first coordinate axis determined bydiffering values of the first parameter and a second coordinate axisdetermined by differing values of the second parameter, and determiningthe values of the parameters from a point on the display, the pointbeing derived from a contour associated with a discrepancy having a verylow value, the value of the first parameter being the coordinate valueof the point along the first coordinate axis and the value of the secondparameter being the coordinate value of the point along the secondcoordinate axis.
 16. The method of claim 15 wherein the first and secondcoordinate axes are respectively n and A in the equation V_(Z)=AZ^(1/n), where Z=depth in the unit and V_(z) is velocity of the waveat depth Z.
 17. The method of claim 15 wherein the method is performedin response to the acoustic wave travelling through a first unit of theearth after having propagated through a second unit of the earth, thefirst and second coordinate axes being respectively V_(o) and k in theequation V_(z) =V₀ +kZ where V_(z) is velocity of the wave at depthZ,V_(o) is velocity at the depth of the first unit closest to the secondunit, and k=the reciprocal of travel time of the wave in the first unit.18. The method of claim 15 wherein the unit includes segments havingdiffering characteristics, each having a differing one of said contours,said differing contours having values overlapping on the display, anddetermining the values of the parameters in one of the segments from apoint on the display where the differing contours overlap.
 19. Themethod of claim 15 wherein the two variables describe a velocity vs.depth characteristic.
 20. The method of claim 15 wherein the twovariables describe a travel time vs. depth characteristic.
 21. Themethod of claim 15 wherein the two variables describe a slowness vs.depth characteristic.
 22. A method of approximately determining firstand second parameters of a theoretical equation having two parametersand two variables describing how an acoustic wave travels as a functionof depth through a unit of the earth, the acoustic wave travellingthrough the unit to a sonic/seismic transducer, comprising responding tosonic/seismic data derived by the transducer responding to the acousticwave by establishing data representing at least one two dimensionalcontour, each contour representing a value of acoustic wave traveldiscrepancy from the theoretical equation of the acoustic wave traveltime to the transducer through the first unit, deriving datacommensurate with the at least one contour as if the contour werepresented on a two coordinate axis display having a first coordinateaxis determined by differing values of the first parameter and a secondcoordinate axis determined by differing values of the second parameter,and determining the values of the parameters from data valuescommensurate with a point on such a display, the point being derivedfrom data corresponding with a contour associated with a discrepancyhaving a very low value, the value of the first parameter beingcommensurate with the coordinate value of the point on such a displayalong the first coordinate axis and the value of the second parameterbeing the coordinate value of the point on such a display along thesecond coordinate axis.
 23. The method of claim 22 wherein the first andsecond coordinate axes on such a display are respectively n and A in theequation V_(z) =AZ^(1/n), where Z=depth in the unit and V_(z) isvelocity of the wave at depth Z.
 24. The method of claim 22 wherein themethod is performed in response to the acoustic wave travelling througha first unit of the earth after having propagated through a second unitof the earth, the first and second coordinate axes on such a displaybeing respectively V_(o) and k in the equation V_(z) =V₀ +kZ, whereV_(z) is velocity of the wave at depth Z,V_(o) is velocity at the depthof the first unit closest to the second unit and k=the reciprocal oftravel time of the wave in the first unit.
 25. The method of claim 22wherein the unit includes segments having differing characteristics,each having a differing one of said contours on such a display, saiddiffering contours having values overlapping on such a display anddetermining the values of the parameters in one of the segments from apoint where the differing contours on such a display overlap.
 26. Themethod of claim 22 wherein the two variables describe a velocity vs.depth characteristic.
 27. The method of claim 22 wherein the twovariables describe a travel time vs. depth characteristic.
 28. Themethod of claim 22 wherein the two variables describe a slowness vs.depth characteristic.
 29. A method of obtaining a good fit between firstand second plots of travel times vs. depth of first and second acousticwaves propagating in a unit of the earth, the acoustic waves travellingfrom the unit to a seismic transducer, the plots being on a twodimensional display having a travel time coordinate and a depthcoordinate, comprising the step of translating the plots relative toeach other along the travel time coordinate so the plots substantiallycoincide on the display.
 30. The method of claim 29 wherein one of theplots has a known velocity vs. depth characteristic, and determining thevelocity vs. depth characteristic of the other plot from the knownvelocity vs. depth characteristic of the one plot when the first andsecond plots substantially coincide on the display.